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HIGHLIGHTS 


•  Reduced  order  model  for  cells  with  phase  change  electrodes. 

•  Two  phase  coexistence  described  using  a  concentration  dependent  diffusivity. 

•  Easy  implementation  of  degradation  to  study  the  cycle  life. 

•  Model  captures  experimentally  observed  signatures  of  phase  transition  and  aging. 
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The  electrochemical  model  has  the  potential  to  provide  a  robust  and  accurate  battery  management 
system,  but  is  not  the  preferred  choice  as  it  involves  solving  non-linear,  coupled  partial  differential 
equations.  In  the  present  work,  a  model  order  reduction  of  the  complete  electrochemical  model  for  a 
lithium  ion  cell  with  phase  change  electrodes  is  carried  out.  The  phase  change  phenomenon  is  described 
using  a  simple,  concentration-dependent  diffusivity  derived  from  mixture  rules.  This  reduced  order 
model  (ROM)  is  validated  with  experimental  data  from  literature.  The  applicability  of  the  model  to 
capture  the  atypical  behavior  of  the  phase  change  electrode  system  is  demonstrated.  Using  the  cell 
response  from  ROM,  charge— discharge  asymmetry  and  path  dependence  in  a  lithium  iron  phosphate 
(LFP)  cell  are  explored  in  detail.  In  addition,  side  reaction  kinetics  and  solid  electrolyte  interphase  for¬ 
mation  are  included  in  the  ROM  framework  to  enhance  its  capability  to  predict  cell  aging.  The  model  is 
used  to  investigate  capacity  losses  occurring  in  a  phase  change  electrode  cell.  Insights  from  these  results 
are  used  to  suggest  cell  operating  guidelines  for  maximizing  utilization. 

©  2014  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

With  the  advent  of  alternate  energy  automobiles  such  as  electric 
and  hybrid  electric  vehicles,  into  the  commercial  market,  there  is  an 
ever  increasing  need  for  a  robust  Battery  Management  System 
(BMS).  Within  the  numerous  functions  of  a  BMS,  the  most  impor¬ 
tant  is  the  accurate  prediction  of  battery  state  and  life.  Conventional 
state  estimators  make  use  of  equivalent  circuit  models  (ECMs)  [1], 
In  spite  of  their  high  speed  and  accuracy,  these  state  estimators  are 
weighed  down  by  their  non-predictive  nature  and  use  of  a  large 
number  of  empirical  parameters.  In  the  quest  for  longer  lasting 
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batteries,  better  electrodes  involving  phase  change  [2]  or  engi¬ 
neered  composite  materials  [3]  are  developed.  As  ECM  based  esti¬ 
mators  work  best  within  the  range  of  data  for  which  they  are  tuned, 
it  is  difficult  to  incorporate  the  atypical  behavior  inherent  to  the 
novel  electrode  materials  using  these  methods. 

Physics  based  models  are  more  generic  than  the  equivalent 
circuit  models  as  the  processes  occurring  within  the  cell  are  rep¬ 
resented  accurately  in  the  former.  In  this  context,  the  electro¬ 
chemical  model  (EM)  proposed  by  Doyle  et  al.  [4  is  applied  to  a 
wide  variety  of  electrochemical  systems,  including  the  lithium  ion 
cells.  The  strength  of  this  approach  lies  in  its  applicability  to  novel 
electrode  materials  that  show  atypical  behaviors,  like  the  multi¬ 
phase  coexistence.  Intercalation/deintercalation  of  lithium  in 
certain  electrode  materials  can  lead  to  the  formation  of  one  or  more 
new  phases.  This  phase  changing  property  of  the  electrode  most 
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often  manifests  as  a  non-monotonic  response  of  the  cell  voltage. 
This  feature  attracts  considerable  research  intrigue  and  hence  is  a 
widely  studied  area.  For  certain  electrodes  such  as  lithium  iron 
phosphate  (LFP)  the  two  phase  coexistence  results  in  a  constant 
discharge  voltage  [2]  which  holds  considerable  industrial  interest. 
Widely  used  negative  electrode  material,  graphite,  is  also  known  to 
transform  to  a  number  of  phases  during  cycling  [5]  adding  to  the 
importance  of  this  phenomenon  in  the  electrochemical  model. 

A  variety  of  approaches  are  adopted  in  the  literature  to  include 
the  phase  change  phenomenon  within  the  electrochemical  model. 
These  range  from  simple  physical  models  which  assume  the  phase 
formation  to  be  isotropic  to  complex  thermodynamics  based 
models  which  take  into  account  the  anisotropy  in  phase  formation. 
In  one  of  the  pioneering  efforts,  Srinivasan  and  Newman  [6] 
represent  the  new  phase  formation  using  the  shrinking  core 
physics.  Instead  of  an  active  material  sphere  with  a  single  phase, 
this  model  assumes  a  shell  and  core  with  the  phase  boundary  at  an 
equilibrium  concentration.  This  approach  is  applied  to  model  the 
discharge  of  an  LFP  half-cell  and  an  LFP/Carbon  (C)  full  cell  [7—9], 
The  extension  of  the  shrinking  core  approach  to  model  the  charge 
part  of  the  cycle  is  proposed  via  a  generalized  moving  boundary 
model  by  Khandelwal  et  al.  [10],  Using  this  approach,  features  of 
the  phase  transformation  in  LFP,  namely,  an  asymmetry  in  the 
charge  and  discharge  cycles,  path  dependence  of  the  response  and 
the  tangential  phase  propagation  within  the  LFP  material  are 
captured.  Motivated  by  the  experimental  evidence  of  anisotropic 
lithium  insertion  [11,12],  phase  field  models  have  been  proposed 
[13—16].  These  models  are  most  useful  to  understand  detailed  as¬ 
pects  of  the  phase  transformation  physics  but  can  prove  difficult  to 
translate  to  onboard  algorithms.  A  major  difficulty  with  the  phase 
change  system  is  to  track  multiple  phase  boundaries  formed  during 
cell  cycling.  This  is  addressed  in  the  generalized  moving  boundary 
approach  [10]  within  the  shrinking  core  family  of  models,  while  it  is 
obtained  as  wave  solution  of  the  phase  field  parameter  in  the  phase 
field  family  [13]  of  models,  both  involving  complicated  solution 
schemes. 

In  an  alternate  approach,  instead  of  considering  a  physical 
boundary,  researchers  make  use  of  a  concentration  dependent 
diffusion  coefficient  [17,18],  The  concentration  dependence  mani¬ 
fests  as  a  thermodynamic  factor/activity  correction  multiplying  the 
constant  diffusion  coefficient.  But,  as  electrolyte  effects  are  not 
considered,  this  model  is  applicable  only  at  rates  until  1 C.  Safari  and 
Delacourt  [19]  also  make  use  of  a  concentration  dependent  diffu¬ 
sion  coefficient  albeit  of  empirical  origin.  The  advantages  of  this 
model  include  a  simplified  approach  to  modeling  the  phase  change 
whilst  considering  the  highly  resistive  nature  of  the  LFP  material 
using  a  resistive-reactant  model  [20],  The  model  is  validated  with 
experiments  on  an  LFP  half-cell  [19]  and  an  LFP/C  commercial  cell 
[21].  Although  the  model  successfully  captures  and  explains  fea¬ 
tures  such  as  charge— discharge  asymmetry  and  path  dependence 
in  LFP,  the  transport  coefficients  are  obtained  by  fitting  the  model 
results  to  experiments  increasing  the  dependency  on  cell  data. 

A  lithium  ion  cell  undergoes  degradation  from  a  variety  of 
sources  [22].  Electrodes  are  known  to  change  in  volume  during  the 
insertion/deinsertion  process  leading  to  fatigue  and/or  cracking 
while  certain  positive  electrode  materials  are  known  to  dissolve 
into  the  electrolyte  [23],  The  most  common  and  well-studied 
mechanism  is  the  occurrence  of  side  reactions  within  the  cell. 
This  involves  reduction  of  the  electrolyte  at  the  negative  electrode 
leading  to  the  formation  of  a  Solid  Electrolyte  Interphase  (SEI).  The 
electrochemical  model  has  been  used  extensively  to  predict  the 
degradation  of  lithium  ion  cells.  Studies  incorporating  just  the 
solvent  decomposition  reaction  [24]  or  combined  with  the  SEI 
formation  [25]  have  been  carried  out  in  the  past.  The  use  of 
empirical  equations  to  model  the  capacity  loss  [26]  or  the  use  of 


first  principles  based  equations  describing  the  SEI  growth  [27]  are 
incorporated  in  the  aging  models.  Models  predicting  calendar  life 
[28]  or  cycle  life  [27]  of  lithium  ion  cells  are  also  studied.  Safari  et  al. 
[25]  present  a  multimodal  aging  model  for  the  lithium  cobalt  oxide 
(LCO)  cell,  valid  for  cycle  life  and  both  constant  voltage  and  OCV 
storage.  In  all  of  the  above  studies  the  degradation  model  is 
included  in  one  form  or  the  other  but  is  valid  for  cells  with  single 
phase  active  material  electrodes.  In  terms  of  modeling  the  degra¬ 
dation/aging  of  a  phase  change  electrochemical  cell,  Safari  and 
Delacourt  [29]  model  the  aging  of  a  commercial  LFP/C  cell.  The 
phase  change  electrode  is  modeled  using  the  earlier  approach  [19] 
of  a  concentration  dependent  diffusivity  and  the  resistive-reactant 
model.  The  degradation  is  included  in  the  form  of  an  SEI  evolution 
equation  while  considering  solvent  diffusion  as  the  limiting  step  to 
the  SEI  formation.  Although  it  incorporates  phase  change,  degra¬ 
dation  and  uses  the  single  particle  model  approach,  the  applica¬ 
bility  of  this  model  is  limited  to  low  and  nominal  rates  as  electrolyte 
effects  are  not  included. 

In  order  to  make  the  electrochemical  model  amenable  for  on¬ 
board  use,  it  needs  to  be  combined  with  a  certain  order  reduction 
methodology.  Also,  cell  aging  studies  involve  long  time  simulations 
of  the  electrochemical  degradation  models.  This  proves  time 
consuming  for  a  high  number  of  cycles.  Therefore,  most  often  re¬ 
searchers  try  to  incorporate  simplified  models  which  can  reduce 
the  simulation  run  time.  The  single  particle  model  (SPM)  proposed 
by  Subramanian  et  al.  [30]  and  applied  by  Safari  and  Delacourt  [19] 
to  study  degradation  of  LFP/C  cells,  is  a  widely  used  approach.  It  is, 
however,  found  to  be  valid  only  until  a  current  rate  of  1C.  A  number 
of  reduced  order  models  have  been  proposed  by  various  groups. 
These  reduced  order  models  include  but  are  not  limited  to  the 
following:  enhanced  single  particle  model  [31,32],  reformulated 
models  [33],  state  variable  models  34]  and  proper  orthogonal 
decomposition  (POD)  based  models  [35j.  Of  these  models,  apart 
from  the  single  particle  model,  the  rest  make  use  of  various 
complicated  mathematical  techniques  for  model  order  reduction. 
The  state  variable  models  need  a  post  processing  step  such  as 
residue  grouping  which  involves  sophisticated  mathematics.  In 
addition,  the  reformulated  model  and  the  POD  based  model  require 
an  initial  guess  of  variables  and  initial  data  for  ensemble,  respec¬ 
tively.  In  the  enhanced  SPM,  the  physical  basis  of  the  SPM  is 
maintained  while  the  electrolyte  effects  are  included  via  parabolic 
[32]  or  Galerkin  approximations  [31], 

In  an  alternate  approach  to  reduce  the  model  order,  Kumar  [36] 
propose  a  completely  physics  based  methodology  developed  using 
volume  averaging.  The  solid  phase  equations  and  the  electrolyte 
profiles  are  obtained  consistently  using  profile  based  approxima¬ 
tions  following  the  volume  averaging.  This  reduced  order  model 
(ROM)  approach  has  a  greater  range  of  applicability  in  terms  of 
current  rates,  when  compared  to  the  SPM  and  a  higher  predict¬ 
ability  than  the  numerical  model  order  reduction  techniques.  On 
the  same  lines,  Kumar  et  al.  [37  propose  an  algebraic  ROM  devel¬ 
oped  for  linearly  varying  currents  motivating  that  any  complicated 
current  profile  can  be  built  of  smaller  linear  current  segments. 
These  studies  are  however,  confined  to  single  phase  electrode  cells. 
Also,  prediction  of  degradation  upon  cycling  or  storage  is  not 
studied  using  this  model. 

A  model  for  a  phase  change  system  should  incorporate  the 
phenomenon  in  a  simple  yet  generic  way,  should  be  able  to  predict 
the  cell  capacity  loss  due  to  aging,  should  not  require  extensive 
tuning  to  experimental  data  and  should  be  amenable  for  onboard 
implementation  in  devices.  Although  the  electrochemical  model  is 
an  ideal  approach  to  predict  accurate  cell  response,  to  incorporate  it 
in  a  BMS  algorithm  a  reduced  order  model  is  required.  The  aim  of 
the  present  work  is  to  develop  a  physics  based  reduced  order 
electrochemical  model  which  can  predict  the  response  of  fresh  and 
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aged  cells  with  phase  change  active  materials.  To  our  knowledge, 
such  a  generic  ROM  incorporating  both  phase  change  and  degra¬ 
dation  has  not  been  carried  out.  The  approach  for  model  order 
reduction  is  based  on  volume  averaging  and  profile  based  ap¬ 
proximations  [36].  In  the  present  work,  a  concentration  dependent 
diffusivity  developed  from  standard  mixture  rules  is  used  which 
does  not  involve  any  empirical  coefficients.  The  side  reaction  cur¬ 
rent  is  described  using  the  Tafel  kinetics  and  is  used  in  the  evolu¬ 
tion  equation  of  the  SEI  film. 

The  manuscript  is  structured  as  follows:  The  problem  descrip¬ 
tion  and  the  model  equations  are  discussed  in  the  following  sec¬ 
tion.  Section  3  consists  of  the  results  and  discussion  in  which  the 
validation  of  the  fresh  as  well  as  aged  cell  is  presented.  Following 
this,  the  ability  of  the  model  to  predict  the  observed  atypical 
behavior  of  the  phase  change  system  is  explored.  Subsequently, 
analysis  of  the  LFP/C  cell  aging  is  presented.  Based  on  the  model 
results,  operating  guidelines  are  proposed  to  minimize  the  lithium 
loss. 

2.  Model 

2.3.  Problem  description 

In  the  present  work  a  reduced  order  model  is  proposed  for  an 
LFP/C  cell.  The  model  order  reduction  is  along  the  lines  of  the 
methodology  proposed  by  Kumar  [36],  A  schematic  of  the  cell  is 
given  in  Fig.  1.  To  simplify  the  model  and  to  differentiate  the  effect 
of  the  phase  change  phenomenon,  only  the  positive  electrode  is 
considered  to  undergo  phase  change.  The  negative  electrode,  under 
the  present  range  of  the  lithiation  studied,  is  assumed  to  have  a 
single  phase.  The  subscripts  n,  p,  1  and  2  are  used  to  denote  the 
negative  electrode,  positive  electrode,  solid  phase  and  the  elec¬ 
trolyte  phase  parameters,  respectively.  The  detailed  nomenclature 
of  the  variables  used  is  provided  in  the  list  of  symbols. 


2.2.  Model  equations 

The  electrochemical  model  (EM)  [4]  is  based  on  the  concen¬ 
trated  solution  theory  for  porous  electrode.  The  equations  used  to 
describe  the  system  consist  of  mass  and  charge  balances  in  both  the 
solid  phase  (electrode)  and  the  electrolyte  phase.  A  list  of  the 
governing  equations  and  the  corresponding  boundary  conditions 
involved  in  EM  are  provided  in  Table  1. 


Table  1 

A  compilation  of  governing  equations  and  the  boundary  conditions  used  in  the 
complete  electrochemical  model. 


Solid  phase  mole  balance 

=  P  s  where  k  =  n,  p 

Boundary  conditions  (BCs) 


=J  k, 


-Dlk^|r=0  =  ° 


cik(r.O)  =  Ciko 

Solid  phase  charge  balance 

If  ^  =  ~ak^Jk 

Negative  electrode  BCs 

_(rln^lrlx=o  =  K*) 
_<Tln^r|x=Ln  =  0 


Positive  electrode  BCs 

_  d<l>^  \ 

^IP  9x  \x=Ln+Ls 
901 P 


:  0 


_<Tip-g?  I x=l  ~ 

Electrolyte  phase  mole  balance 

e2k^§r  =  If  +  fln0  -  f+l/k 

Negative  electrode  BCs  Positive  electrode  BCs 

-^2n^  |x=0  =  0  -^2p-|r  lx=L„+Ls  =  ^2ip 

~D2n^  lx=L„  =  ^2in  ~D2p^§r  |X=L  =  ® 

c2i\\x=Ld  =  C2in  C2plx=tn+Ls  =  C2ip 


Total  charge  balance  (solid  +  electrolyte) 

.-^-K,k^  +  2K2n¥(l-t+)a|A  =  / 


-°lk  i 
Negative  electrode  BCs 

~K2nT5xL  lx=0  =  0 
-k2n%|x=lo  =/(t) 
<f>2n\x=ia  =  02in(O 


Positive  electrode  BCs 

=L.+L,='W 

=L  z 


a^2p  I 

-K2psr  I; 

-Ko  ^1 

K2p  ax  I; 


:  0 


^2plx=Ln+ts  —  02ip(O 


Separator  BCs 

|x=Ln  =  Qhn 

-02s ^  lx=Z.n+Ls  =  ^2ip 
C2s|x=  Ln  =  C2in 
C2slx=Ln+Ls  =C2iP 


Separator  BCs 

-K2S^U  =Kt) 

~K2s  %■  x=Ln+rs  = 
<hs\x=ln  =02in(O 

feL=Ln+Ls  =fep(t) 


The  governing  equations  are  partial  differential  equations 
(PDEs),  wherein,  the  dependent  variable  has  both  spatial  and 
temporal  dependence.  To  reduce  the  model  order,  in  the  present 
work,  the  philosophy  followed  is  to  overcome  the  spatial  depen¬ 
dence  [36],  The  simplest  procedure  to  overcome  spatial  depen¬ 
dence  of  variables  is  to  work  with  their  average  values.  Although, 
this  directly  reduces  the  PDEs  to  ordinary  differential  equations  this 
step  also  leads  to  loss  of  spatial  information  of  the  system  variable. 
Therefore,  the  volume  averaged  variables  are  coupled  with  suitable 
profile  approximations  to  gain  back  the  spatial  information.  The 
steps  followed  are  as  given  below: 

To  begin  with,  the  average  of  any  system  variable  f[x,t)  across  a 
region  of  thickness  L  with  a  surface  area  S,  in  the  Cartesian  co¬ 
ordinates  is  defined  as, 

L 

(f(t)>  =  ^I//(^t)Sdx!  (1) 

o 

where  (  )  denotes  the  volume  averaged  quantity.  The  average  of  a 
variable  in  the  radial  coordinates  (denoted  by  an  overbar)  for  a 
particle  of  radius  R  is  defined  as, 

R 

7<t)-dn5  /4"r2/(r,dr-  <2> 

3  o 

The  primary  objective  of  the  model  is  to  estimate  the  cell  voltage 
which  is  the  potential  drop  across  the  positive  and  the  negative 
electrodes. 


plx— t  01nlx=O 


Fig.  1.  A  schematic  representation  of  the  cell. 


Kzell  —  01 


(3) 
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The  solid  phase  potentials  (</>ip,  </>ln)  can  be  calculated  from  the 
inverted  Butler— Volmer  equations, 


that  the  overpotential  driving  the  side  reaction  is  affected  by  the 
Ohmic  drop  across  this  film.  The  film  resistance  is  calculated  as, 


^iplx=t  ~  02p|x=L  +  up(c*p)  +— ^sinh  1  ^  2jp0L) 

^ln|x=o  =  02n|x=O  +  ^n(Csn)  H - sinh  1  Jj1  (5) 

where,  <p2p\x=L  &  02n|x=o  are  the  electrolyte  phase  potentials  within 
the  positive  and  the  negative  electrodes,  at  the  respective  current 
collector  ends  and  jp;  &  jm-  are  the  corresponding  rates  of  the 
intercalation  reactions.  Up  and  Un  are  the  open  circuit  potentials  of 
the  positive  and  the  negative  electrode,  respectively,  and  are 
functions  of  the  lithium  surface  concentrations  in  the  electrode 
particles  (csp,  csn).  The  surface  concentrations  are  also  coupled  to 
the  electrolyte  phase  through  the  reaction  rate  prefactors  (Jn o,  jpo) 
given  by, 


JnO  —  kn(Clmaxn  -  Csn)  Cjffcjn,  (6) 

jpo  =  k p(ctmaxp  -  csp)°'5c^c^.  (7) 


kn,  kp  are  the  reaction  rate  constants,  cimaXn,  cimaxp  are  the 
maximum  concentrations  within  the  solid  phase  and  C2n,  C2P  are  the 
electrolyte  phase  concentrations,  in  the  negative  and  positive 
electrodes  respectively. 

In  the  present  work,  the  side  reaction  involving  electrolyte 
reduction  at  the  negative  electrode  is  considered  to  be  the  main 
source  of  degradation.  The  total  current  at  the  negative  electrode  is 
split  into  the  component  responsible  for  the  lithium  intercalation/ 
deintercalation  reaction  and  the  one  responsible  for  the  side  re¬ 
action.  This  is  explained  using  a  schematic  in  Fig.  2. 


=>jn(x,  t)  =jni(x,  t)  +js{x,  t)  (8) 

The  side  reaction  rate  can  be  modeled  using  Tafel  kinetics  [24] 
as, 

Js(x,  t)  =  ^jjexp(-  ^  (</>ln  -  <t>2 n  -  Us -jn(x,  (9) 

The  side  reaction  results  in  the  formation  of  a  resistive  film  of 
reaction  products  on  the  electrode  surface  (SEI).  It  should  be  noted 


Rf  =  A 


A 

Kf  ’ 


(10) 


where,  Kf  is  the  conductivity  of  the  SEI  film  while  <5f  is  the  film 
thickness  which  is  obtained  from, 


d<5f  _  jsMf 


(H) 


That  is,  evolution  of  the  film  depends  on  its  molecular  weight 
(Mf),  density  (pf)  and  rate  of  the  side  reaction.  To  take  the  degra¬ 
dation  into  account,  Eq.  (5)  is  modified  to, 


01n|x=O  -  ^2n|x=o  +  Un(CSn)  +  Jn|x=ofRf  lx=0 

+  ^Sinh-1  pn|x=o-jS|x=0^  (12) 

F  \  2 ]n0  J 

In  the  absence  of  side  reaction  at  the  positive  electrode,  Eq.  (4) 
remains  unchanged.  In  this  case,  the  intercalation  reaction  cur¬ 
rent  Opi(x,t))  is  equal  to  the  total  current  (/p(x,t)),  that  is, 

^tplx— t  =  <feP|x=i  +  UpM  +  — ^-sinfT1  •  (13) 

To  determine  the  cell  potential  using  Eqs.  (12)  and  (13),  the 
unknown  variables  are  the  active  material  surface  concentration  of 
lithium  and  the  electrolyte  phase  potentials  at  the  current  collec¬ 
tors.  The  procedure  to  obtain  these  variables  is  outlined  next. 


2.2.1.  Obtaining  the  active  material  surface  concentrations 

As  a  first  step,  the  volume  averaging  procedure  given  by  Eq.  (1) 
is  applied  to  the  solid  phase  charge  balance  Eq.  (refer  Table  1)  for 
the  negative  electrode  region.  Upon  using  the  boundary  conditions, 
this  equation  reduces  to 

I(t)  =  anLnF(jn)(t),  (14) 


where  (jn>  (t)  is  the  volume  averaged  reaction  flux  and  1(f)  is  the 
applied  current.  Extending  this  to  the  positive  electrode  and  rear¬ 
ranging  the  equations  gives  the  relation  between  the  average  re¬ 
action  rates  and  the  applied  current. 


0'n>(t) 


Jit )_ 

Gnbn  F 


(15) 


From 

positive 

electrode 


Fig.  2.  A  schematic  representation  of  the  SEI  formation  on  an  active  material  particle 
surface. 


WW’W  <161 

The  surface  concentration  of  lithium  in  the  active  material 
particles  is  obtained  from  the  solid  phase  mass  balance  equation. 

In  the  present  work,  the  positive  electrode  active  material  is 
considered  to  undergo  phase  change  during  the  insertion/dein¬ 
sertion  of  lithium.  The  key  step  towards  modeling  a  phase  change 
system  is  to  describe  the  solid  phase  diffusion  accurately.  During  a 
major  portion  of  the  charge/discharge  cycle  of  the  cell,  multi  phases 
can  coexist.  The  diffusion  within  the  active  material  is  expected  to  be 
affected  due  to  the  presence  of  this  coexistence.  In  the  past,  the 
diffusion  process  has  been  described  using  a  shrinking-core 
approach  [2,6],  A  couple  of  experimental  studies  [11,12]  of  the 
lithium  insertion/deinsertion  in  LFP  materials  provide  evidence  that 
lithium  diffusion  occurs  preferentially  along  the  b  direction  in  the 
crystal  plane  and  the  phase  boundary  propagates  anisotropically 
throughout  the  particle.  Therefore,  solid  phase  diffusion  models  with 
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anisotropic  diffusivity  [17—19]  are  studied  in  lieu  of  the  isotropic 
shrinking  core  model.  In  the  present  work,  this  phenomenon  is 
described  using  a  concentration  dependent  diffusion  coefficient 
developed  assuming  the  system  to  be  a  mixture  of  two  phases. 

In  this  regard,  the  arguments  used  in  the  coupling  of  the  pseudo 
2D  model  are  applied  to  this  case.  In  the  pseudo  2D  model  for  a 
porous  electrode,  it  is  assumed  that  the  electrode  is  a  homogeneous 
mixture  of  both  the  solid  and  electrolyte  phases.  This  enables  the 
model  to  couple  the  solid  phase  concentration,  obtained  by  solving 
the  solid  mass  balance  equations  (in  radial  coordinates),  to  the 
electrolyte  equations,  at  any  point  along  the  thickness  of  the  cell.  In  a 
phase  change  system,  irrespective  of  the  mode  of  the  phase 
boundary  propagation,  the  active  material  particle  consists  of  the 
two  coexisting  phases.  A  schematic  of  the  arguments  presented  here 
is  provided  in  Fig.  3.  From  the  schematic  it  can  be  seen  that,  along  the 
radius  of  the  particle,  one  encounters  a  mixture  of  both  the  phases. 

Accordingly,  the  diffusivity  of  the  phase  change  active  material 
particle  can  be  expressed  as, 

^IPmix  =  XP^lpa  +  XpP^lpP  (17) 

where,  xpo,,  xpp  are  the  mole  fractions  of  the  lithium-deficient  phase, 
a  and  the  lithium  rich  phase,  0,  respectively  and  Dipa,  Dipp  are  their 
diffusion  coefficients.  Xipa  and  Xipp  are  given  by  the  following 
equations, 


cipP  (cip) 
^-lpP  —  clpix 


(18) 


xipp  —  1  xipct  (19) 

It  can  be  noted  that  the  diffusivity  varies  with  the  average 
concentration,  cip  while  taking  into  account  the  equilibrium  con¬ 
centrations  and  diffusivities  of  the  individual  phases  (a  and  0) 
present  in  the  system.  When  cip  is  equal  to  either  cipa  or  cipp,  the 
diffusion  coefficient  is  equal  to  the  pure  phase  diffusivity,  D^a  and 
Dip,  respectively. 

The  solid  phase  mass  balance  for  the  positive  electrode  under¬ 
going  phase  change  can  be  written  as 


^!p  =  JLA(r2Dlp  ^ 

at  r2  dr  \  lp™-  dr  ) 


(20) 


o 

Isotropic  phase 
change 

» 
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Fig.  3.  A  schematic  representation  of  the  isotropic  shrinking  core  approach  and  the 
anisotropic  phase  propagation  approach  adopted  for  a  phase  change  system  in  liter¬ 
ature.  Irrespective  of  the  model  of  phase  formation,  in  the  present  work,  the  system  is 
considered  to  be  equivalent  to  a  mixture  of  two  phases. 


With  the  BCs:  -Dlp  =joi,  -Dlp 

iiw  ar  |r=R  J P”  IPn 


*l)Flr=0 


=  0 


(21) 


while  that  for  the  negative  electrode  undergoing  degradation  is 
given  by 


9cln  19  /  2n  ®cln 


9 1  r2  dr 


dr 


(22) 


With  the  BCs:  -  Dln^  |r=Rn  =  jni,  -Dln9^|r=0  =  0  (23) 

Applying  the  volume  averaging  [30]  given  by  Eq.  (2)  to  the  solid 
phase  mole  balance  Eqs.  (20)  and  (22)  and  using  the  boundary 
conditions,  leads  to  the  equation, 


(24) 


It  can  be  noted  that  the  above  equation  gives  the  concentration 
of  lithium  averaged  over  one  single  particle.  This  equation  is, 
therefore,  subjected  to  averaging  across  the  electrode  volume  Eq. 
(1)  to  obtain, 


d(clk>  3  (/  ki  > 


dt 


Ru 


(25) 


The  average  reaction  rate  (jki)  is  equal  to  (jp)(t)  for  the  positive 
electrode  and  can  be  replaced  by  the  relation  to  the  applied  current 
given  by  Eq.  (16).  For  the  negative  electrode,  the  relation  given  by 
Eq.  (8)  is  used  to  obtain  the  relation  with  the  total  current. 

To  obtain  the  concentration  gradient  within  the  active  material 
particle,  the  quartic  profile  described  by  Subramanian  et  al.  [30] 
and  Kumar  [36]  is  found  to  be  most  suitable.  Following  the  pro¬ 
cedure  outlined  therein,  the  radial  gradient  of  the  lithium  within 
the  active  materials  can  be  obtained  by  solving, 


dc 


lkr 


30DikCikr  45jki 


dt 


Ri 


2  Rl 


(26) 


wherein,  D ik  =  D1Pmx  for  the  positive  electrode  and  is  given  by  Eq. 
( 17).  These  equations  upon  solving  can  be  used  to  obtain  the  surface 
concentrations  in  both  the  electrodes 


Csk 


-p  R\dki  ,  8RkClkr 

clk  ~  -TFr5 - F  - 


35Di 


35 


(27) 


2.2.2.  Obtaining  the  electrolyte  phase  potentials 

To  obtain  the  electrolyte  phase  potentials  at  the  current  col¬ 
lector  ends,  the  total  charge  balance  listed  in  Table  1  is  used.  This 
equation  includes  contribution  from  the  solid  phase  charge  balance 
and  is  coupled  with  the  electrolyte  phase  concentration. 

To  take  into  account  the  spatial  variation  of  the  solid  phase 
potential  (<jqk).  the  local  reaction  rate  is  assumed  to  be  equal  to  its 
volume  average  value. 

J'k(x,t)  ~  Ok)(t)(k  =  n,p),  (28) 

Substituting  the  value  of  average  reaction  rate  instead  of  the 
local  in  the  solid  phase  charge  balance,  integrating  and  applying  the 
electrode— electrolyte  interface  boundary  condition  leads  to  the 
expression  for  the  electric  field  in  the  solid  phase  of  the  electrodes. 


°Tn 


d<ftin 

9x 


(Ln~ 


*) 


(29) 
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“<riP^?  =  T7(x“(Ln  +  Ls))  (30) 

The  concentration  within  the  electrolyte  is  obtained  by  solving 
the  electrolyte  phase  mass  balance  equations  in  the  electrodes  and 
separator.  The  volume  averaged  forms  of  the  equations  are  ob¬ 
tained  by  using  Eq.  (1 )  and  the  flux  boundary  conditions  to  be, 

=  -<?2in  +  0  “  t+)anLn(jn )  (31 ) 

Lpt'ip — =  <?2ip  +  0  —  f+)ap^-pOp)  (32) 


differential  equations  and  algebraic  equations.  These  equations  are 
incorporated  in  Simulink  through  a  Matlab  function  block.  The 
differential  equations  are  solved  using  integrator  blocks  which  uses 
the  ODE15s  solver.  A  discussion  of  the  results  is  provided  in  the  next 
section. 

3.  Results  and  discussion 

In  the  present  section,  the  results  from  ROM,  which  includes 
multiphase  behavior  of  electrode  materials  and  a  kinetic  SEI  growth 
model  to  study  the  cycle  life  of  cells,  are  presented.  The  first  part 
deals  with  the  discussion  of  fresh  cell  results  followed  by  the  cell 
aging  studies. 


V2S^  =  q2in-q2ip  (33) 

In  order  to  obtain  the  concentration  profiles  within  the  elec¬ 
trolyte  phase,  the  local  value  of  the  gradient  of  diffusive  concen¬ 
tration  flux  (within  the  electrolyte)  is  substituted  with  its  average 
value. 


In  Eq.  (34),  k  takes  the  values  n  and  p  for  the  negative  and  the 
positive  electrode,  respectively.  Integrating  these  equations  in  both 
the  electrodes  and  the  separator  regions  whilst  using  the  concen¬ 
tration  boundary  conditions,  the  electrolyte  phase  concentrations 
are  obtained  as 


C2n(*,  t) 


C2in  (t)  + 


Q  2in(t) 

2LnD2n 


(35) 


c2p(x,  t)  =  c2ip(t) 


*?2ip ( 0  / j2 

2LpD2p  \  P 


(36) 


C2s(*,t) 


C2in(t)  - 


<?2in(t) 
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D2P 

{ Q2in(t)  ~  QlipY) 


(X  -  Ln) 


LD 


2s 


(X-LdY 


(37) 


where,  C2m,  C2iP,  tfein  and  q2iP  are  concentration  and  flux  variables  at 
the  electrode— electrolyte  interface.  They  are  determined  by  sum¬ 
mation  of  the  total  concentration  in  the  electrolyte  phase  and  can 
be  obtained  from  Eqs.  (83),  (84),  (78)  and  (79)  in  Kumar  [36]. 

Knowing  the  concentration  profiles  (Eqs.  (35)— (37))  and  the 
relation  between  the  electric  field  in  the  solid  phase  and  the 
applied  current  (Eqs.  (29)  and  (30)),  the  total  charge  balance  is 
integrated  to  obtain  the  electrolyte  phase  potentials.  At  the  current 
collector  ends  the  values  are, 


^nlx^o  —  ^2in  + 


2RZT  . 

-#-(1  -4)log 


c2n(x  =  0)\  ILn 


c2in 
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2k- 


2n 


(38) 


i  2RgT  „  N1 

02p I X=L  ~  </)2ip  H  p  (1  — 


c2p  (x  =  L) 

c2ip 
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ILn 


2  K 


2p 


(39) 


Using  Eqs.  (27),  (38)  and  (39)  in  Eqs.  (6),  (7),  (12)  and  (13),  the 

solid  phase  potentials  at  the  current  collector  ends  are  obtained. 
These  are  used  to  obtain  the  cell  response  using  Eq.  (3). 

To  summarize  the  discussion  of  the  model  equations,  the  com¬ 
plete  electrochemical  model  for  a  cell  with  a  phase  change  positive 
electrode,  with  degradation,  is  reduced  to  solving  ordinary 


3.1.  Model  validation 

3.1.1.  Reduced  order  model  for  phase  change  system  (js  ~  0) 

To  demonstrate  its  validity,  ROM  for  fresh  cells  (first  cycle)  is 
compared  to  experimental  data  from  literature.  Safari  and  Dela- 
court  [21  ]  report  galvanostatic  charge— discharge  curves  at 
different  rates  for  an  LFP/C  cell.  For  the  present  work,  data  points 
for  two  rates  namely,  C/10  and  1 C  are  used  for  model  validation.  The 
values  of  the  cell  parameters  used  are  same  as  those  reported  by  Ye 
et  al.  [38],  Safari  and  Delacourt  [21]  and  are  listed  in  Table  2.  A 
comparison  of  model  results  with  experimental  data  is  shown  in 
Fig.  4.  Good  agreement  between  the  experimental  and  simulated 
values  is  observed,  especially  at  the  beginning  and  end  of  the 
charge/discharge  cycles.  A  relative  mismatch  is  observed  in  the  mid 
portion  of  the  curves.  This  mismatch  is  seen  to  be  higher  during 
discharge  than  during  the  charge  cycle.  It  is  also  found  to  be 
increasing  with  rate  but  the  extent  of  mismatch  is  found  to  be  less 
than  1.8%.  This  substantiates  the  utility  of  this  phase  change 
reduced  order  model  at  both  low  and  nominal  rates. 

3.1.2.  Reduced  order  model  for  phase  change  system  with 
degradation 

To  study  the  degradation  of  LFP/C  cell,  the  side  reaction  rate  is 
assumed  to  have  a  finite  value  during  the  charging  cycle.  Safari  and 
Delacourt  [39]  report  experimental  cell  response  data  for  various 
cycles  for  an  LFP/C  cell.  As  the  initial  states  (SOCo  =  Cio/cimax)  for 
these  set  of  experiments  are  different,  the  values  of  SOC„o  =  0.79, 
SOCpo  =  0.03  are  used  for  validating  the  cycling  results.  A  com¬ 
parison  of  model  results  with  experimental  data  for  the  first  and 


Table  2 

A  list  of  the  parameters  used  in  the  model. 


S.  no 

Parameter 

LiFeP04 

LiC6 

Separator 

Source 

1 

L  (m) 

70  x  10" 

6  34  x  10-6 

30  x  10-6 

[38] 

2 

r(m) 

36.5  x  10-9  3.5  x  10-6 

[38] 

3 

£2 

0.332 

0.33 

0.54 

4 

e\ 

0.43 

0.55 

5 

Dt  (m2  s-1) 

Eq.  (17) 

3.9  X  to-14 

1381 

6 

D2  (m2  s-1) 

2.8  X  10-10 

1381 

7 

cimax  (mol  nrr3) 

22,806 

31,370 

[381 

8 

Cio/Cl  max 

0.035 

0.82 

9 

C20  (mol  m-3) 

1200 

9 

k  (m  s_1) 

1.4  x  10 

12  8.19  x  10-12 

[381 

10 

k  (S  nr1) 

1.3634 

1.3459 

1.3566 

[381 

11 

S(m2) 

0.1694 

0.1755 

S.  no 

Parameter 

SEI  film 

Source 

1 

js0  (mol  nr3 

s-1) 

3.4  x  10-”  x  0.05 

x  4.541 

[29J 

2 

Pt  (kg  nr3) 

1690 

[291 

3 

Mf  (kg  mol-1 

) 

0.162 

[291 

4 

k{  (S  m-1) 

3.8  X  10~7 

[29] 

5 

14  (V) 

0.38 

[29] 
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Fig.  4.  The  model  results  are  validated  with  the  experimental  data  reported  in  Safari  Fig.  6.  Charge  and  discharge  curves  at  different  rates  are  plotted  for  the  LFP/C  cell 
and  Delacourt  [21]  at  different  rates  for  both  charge  and  discharge  cycles.  under  study.  An  asymmetry  in  the  end  capacities  is  observed. 


the  310th  cycle  is  given  in  Fig.  5.  An  excellent  match  between  the 
model  and  experimental  results  with  <1%  error  is  observed.  This 
benchmarks  the  ROM  approach  used  in  the  present  work  to  model 
degradation  of  a  phase  change  cell  and  demonstrates  the  applica¬ 
bility  of  the  model. 

3.2.  Charge— discharge  asymmetry 

Due  to  multi-phase  coexistence,  an  interesting  feature  of  LFP/C 
cells  is  the  charge— discharge  asymmetry.  In  Fig.  4,  an  asymmetry  in 
the  end  capacities  (utilization)  during  charge  and  corresponding 
discharge  cycle,  is  observed.  To  explore  this  observation  further, 
rate  capability  results  for  the  fresh  LFP/C  cell  are  studied.  Assuming 
a  completely  charged  initial  state,  the  cell  discharge  at  four 
different  rates  is  shown  (solid  lines)  in  Fig.  6.  At  the  same  rates,  the 
cell  response  during  charge  from  a  completely  discharged  state  is 
shown  as  dashed  lines.  At  low  rates,  say  C/10,  the  end  capacity  for 
discharge  and  charge  are  the  same.  While  the  end  capacities  do  not 
vary  significantly  during  discharge,  the  decrease  in  the  end  ca¬ 
pacities  during  charge  is  quite  significant.  This  gives  rise  to  an 
asymmetry  in  the  charge  and  discharge  profiles.  This  asymmetric 
behavior  is  characteristic  of  the  phase  change  electrode  and  has 
been  reported  previously  [10,19].  Thus  it  is  observed  that  the  pro¬ 
posed  model  can  be  used  to  characterize  the  observed  atypical 
features  of  phase  change  cells. 

In  Fig.  7  the  end  capacities  obtained  during  charge  and  discharge 
are  plotted  as  a  function  of  the  current  rate,  it  is  observed  that  at 
any  particular  rate,  the  end  capacity  obtained  via  charging  is  lower 


Fig.  5.  The  model  results  are  validated  with  the  cycle  life  experimental  data  reported 
in  Safari  and  Delacourt  [39]  at  different  discharge  cycles. 


than  that  obtained  via  discharging.  The  difference  in  the  end  ca¬ 
pacities  is  found  to  increase  with  increasing  rate.  The  feature  is 
observed  in  experiments  by  Safari  and  Delacourt  [21], 

To  understand  the  asymmetry  further,  concentration  profiles 
during  charge  and  discharge,  within  the  positive  electrode,  are 
investigated.  The  evolution  of  the  average  concentration  for  a  1C 
discharge  is  shown  in  Fig.  8  as  a  solid  line,  while  that  during  charge 
is  given  by  a  dashed  line.  It  is  observed  that  the  change  in  the 
concentration  during  charge  and  during  discharge  is  the  same,  but, 
the  time  taken  to  do  so,  and  hence  the  capacity,  is  different.  The 
asymmetric  response  is  attributed  to  the  different  values  of  trans¬ 
port  coefficients  between  charge  and  discharge  [40],  In  the 
shrinking  core  approach  [6,10]  this  is  taken  into  account  by 
considering  different  diffusivities  of  the  a  and  /?  phases.  The  present 
model  accounts  for  the  same  behavior  by  the  concentration 
dependent  diffusion  coefficient  (Eq.  (17)). 

3.3.  Path  dependence 

The  heterogeneous  environment  in  the  active  material  particle, 
due  to  the  coexistence  of  two  phases,  affects  the  transport  of 
lithium  during  insertion  differently  than  during  deinsertion.  For 
example,  for  LFP  during  discharge,  a  lithium  rich  phase  (/?)  forms  at 
the  cost  of  the  initial  lithium  poor  phase  (a)  and  the  reverse  occurs 
during  the  charge  cycle.  As  the  diffusion  coefficients  are  different 
for  these  two  phases  [40],  lithium  transport  correspondingly  takes 
shorter  or  longer  depending  on  the  phase  it  travels  through. 

A  complete  and  accurate  prediction  of  path  dependence  of  a  cell 
is  an  important  requirement  for  building  the  BMS  [40],  This  feature 


Fig.  7.  The  end  capacities  of  the  charge  and  discharge  curves  in  Fig.  6  are  plotted  as  a 
function  of  the  current  rates. 
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Capacity  (Ah) 

Fig.  8.  The  average  concentration  of  the  two  phase  system  varying  during  discharge 
and  charge.  The  mismatch  in  the  end  concentrations  is  indicated  by  the  arrow. 


has  earlier  been  explained  using  the  shrinking  core  model  40]  and 
the  generalized  moving  boundary  model  [10],  These  models  have 
phase  boundaries  defined  which  specify  the  diffusion  lengths 
involved  in  different  paths.  In  the  present  work,  in  order  to  make 
the  model  amenable  for  multi-cycling  studies,  no  explicit  phase 
boundary  is  considered.  It  is  therefore  interesting  to  study  if  the 
path  dependence  can  be  captured  by  the  present  description  of  the 
phase  change  system. 

To  demonstrate  the  model's  capability  in  predicting  the  path 
dependent  behavior  of  the  phase  change  cell  (LFP/C),  two  different 
protocols  are  used.  In  the  first  protocol,  the  time  for  the  different 
paths  is  fixed  allowing  the  concentration  to  vary  while  in  the  sec¬ 
ond  protocol,  the  concentration  to  be  attained  at  the  end  of  the  path 
is  fixed. 

3.3.1.  Fixed  time  protocol 

In  the  fixed  time  protocol  to  study  path  dependence,  the 
response  of  the  cell  subjected  to  two  different  paths  is  compared.  In 
the  first  path  (PI ),  the  cell  is  discharged  from  a  completely  charged 
state  for  2000  s.  This  forms  the  first  stage  of  the  path.  In  the  second 
stage  it  is  let  to  rest  for  2  h.  This  long  rest  period  ensures  that  the  cell 
is  equilibrated.  In  the  third  stage,  post  resting,  the  cell  is  discharged 
until  it  reaches  the  cut-off  voltage  of  2  V.  The  second  path  (P2 )  differs 
from  PI  only  in  the  first  stage,  wherein,  the  cell  is  charged  from  a 
completely  discharged  state  for  2000  s.  The  following  two  stages, 
the  rest  and  discharge  steps,  remain  same  as  those  in  PI.  The  details 
of  the  paths  are  summarized  in  Table  3.  The  cell  potential  for  the 
complete  path  is  plotted  as  a  function  of  the  cell  capacity  in  Fig.  9. 
Irrespective  of  the  initial  stages,  in  both  the  cases,  the  cell  is  rested 
for  two  hours,  yet,  a  clear  difference  in  the  end  capacities  is 
observed.  The  discharge  curve  of  P2,  with  an  initial  charging  stage,  is 
found  to  have  a  higher  end  capacity  than  the  other. 

To  gain  a  physical  understanding  of  the  path  dependence,  the 
average  concentration  of  lithium  in  the  active  material  particles,  for 


Table  3 

A  summary  of  the  various  paths  used  in  the  present  work,  to  study  the  path 
dependence  in  LFP/C  cells. 


Path 

Stage  1 

Stage  2 

Stage  3 

Result 

Fixed  time  protocol 

PI 

1C,  t=  2000  s 

0C,  t  =  7200  s 

ic,  Vcut„ff=2V 

Lower  utilization 

P2 

-1C,  t=  2000  s 

0C,  t  =  7200  s 

1C,  Vcut„ff=2V 

Fixed  concentration  protocol 

PI 

1C,  SOC  =  0.5 

0C,  t  =  7200  s 

1C,  VCutoff=2V 

P2 

-1C,  SOC  =  0.5 

0C,  t  =  7200  s 

1C,  VCutoff  =  2V 

Lower  utilization 

P3 

-1C,  SOC  =  0.5 

0C,  t  =  7200  s 

-1C,  Vcutoff=3.5V 

P4 

1C,  SOC  =  0.5 

0C,  t  =  7200  s 

-1C,  Vcut„ff=3.5V 

Lower  utilization 

Fig.  9.  The  simulated  response  of  the  LFP/C  cell  to  two  different  current  paths. 


both  the  paths,  is  plotted  as  a  function  of  the  cell  capacity  in  Fig.  10. 
It  can  be  seen  from  the  figure,  at  the  beginning  of  the  stage  3,  that 
the  average  concentration  in  the  cell  while  undergoing  PI  is  higher 
than  when  undergoing  P2.  This  mismatch  is  the  result  of  the 
asymmetry  in  the  charge  discharge  utilization  as  discussed  in  the 
Section  3.2.  This  in  turn,  leads  to  an  asymmetry  in  the  end  capac¬ 
ities  originating  from  the  initial  paths  taken. 

3.3.2.  Fixed  concentration  protocol 

In  this  protocol,  the  cell  is  subjected  to  four  different  paths.  As  in 
the  earlier  protocol,  all  the  paths  involve  three  stages.  In  the  first 
path  (PI),  the  first  stage  involves  discharging  the  cell  to  a  state  of 
charge  (where  SOC  =  Csp/cimaxp)  of  0.5.  That  is,  instead  of  fixing  the 
time  as  done  in  section  (3.3.1 ),  the  concentration  to  which  the  cell  is 
discharged  is  fixed.  In  the  second  stage,  the  cell  is  rested  for  2  h 
followed  by  a  third  stage,  wherein,  the  cell  is  discharged  from  this 
state  to  a  cut  off  voltage  of  2  V.  The  second  path  (P2),  differs  from 
the  first  path  only  in  the  first  stage.  In  P2,  the  cell  is  charged  from  a 
completely  discharged  state  to  SOC  =  0.5  followed  by  a  2  h  rest  and 
1C  discharge  to  2  V.  Path  3  (P3)  is  the  reverse  of  PI.  In  the  first  stage 
of  P3,  the  cell  is  charged  from  a  completely  discharged  state  to 
SOC  =  0.5.  This  is  followed  by  rest  for  2  h  and  1C  charging  to  a  cut¬ 
off  voltage  of  3.5  V.  The  fourth  path  (P4)  differs  from  P3  only  in  the 
first  stage.  It  involves  a  discharge  to  the  mid  SOC  followed  by  rest 
and  then  charge  to  the  cut-off  of  3.5  V.  The  details  of  the  paths  are 
summarized  in  Table  3.  The  response  of  the  cell  to  all  the  four  paths 
is  plotted  in  Fig.  11.  As  the  mid  SOC  is  reached  at  different  times 
during  charging  and  discharging,  the  cell  response  ends  at  different 
times.  To  exclude  this  effect,  the  time  after  rest,  that  is,  the  begin¬ 
ning  of  the  stage  3  for  all  the  four  paths  is  considered  to  be  initial 
time  (f  =  0)  in  the  inset  of  the  same  figure. 


Capacity  (Ah) 

Fig.  10.  The  variation  of  the  average  concentration  of  the  lithium  in  the  positive 
electrode  during  the  paths  1  and  2. 
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From  the  inset  in  Fig.  11,  a  capacity  loss  in  both  the  discharge  and 
the  charge  cycles  can  be  observed.  Within  the  discharge  curves  (PI 
and  P2),  the  curve  which  has  a  history  of  a  charge  cycle  (P2)  shows  a 
lower  utilization  while  within  the  charge  curves  (P3  and  P4),  the  path 
(P4)  with  a  history  of  a  discharge  cycle  shows  a  lower  utilization. 

To  understand  the  origin  of  this  path  dependence,  the  surface 
and  average  concentrations  in  the  positive  electrode  particles  are 
plotted  as  a  function  of  the  cell  capacity  for  the  two  discharge  paths 
PI  and  P2  in  Fig.  12.  From  the  figure,  it  can  be  noticed  that  at  the  end 
of  the  stage  1,  the  surface  concentrations  reach  0.5  x  clmaXp  but 
upon  resting  in  stage  2  they  relax  to  different  average  concentra¬ 
tions.  These  average  concentrations  form  the  initial  state  for  the 
stage  3.  As  the  stage  3  in  the  PI  begins  with  a  lower  concentration, 
the  corresponding  diffusion  coefficient  is  higher  resulting  in  better 
lithium  transport  and  higher  utilization. 

Above  discussion  leads  to  following  observations:  Partial 
charging  of  phase  change  cells  from  a  completely  discharged  state, 
before  use  (discharge),  can  lead  to  lower  utilization.  On  the  same 
lines  paths  3  and  4  suggest  that  while  charging,  higher  capacity  can 
be  attained  if  the  preceding  stage  is  charging.  That  is,  lower  utiliza¬ 
tion  can  occur  when  charging  a  partially  discharged  cell.  The  analysis 
of  the  model  results  show  that  the  utilization  can  be  maximized  if 
similar  stages  (charge  or  discharge)  succeed  one  another.  The  in¬ 
sights  from  this  path  dependence  study  can  be  used  as  operating 
guidelines  for  cells  such  as  LFP/C  that  exhibit  phase  change. 

3.4.  Cycle  life 

In  this  section,  cycle  life  of  the  LFP/C  cell  is  studied.  As  discussed 
in  Section  2,  the  cell  degrades  due  to  side  reactions  occurring  on  the 


Fig.  12.  The  surface  and  average  concentrations  within  the  positive  electrode  of  the 
LFP/C  cell  varying  due  to  two  different  current  paths. 


negative  electrode  surface  during  aging.  The  reaction  consumes 
lithium  ions  causing  capacity  loss  while  the  products  form  an  SEI 
film  on  the  active  material  particle  surface  increasing  the  cell 
resistance  while  in  operation.  The  proposed  model  has  been  vali¬ 
dated  with  experimental  data  from  the  literature  in  Section  3.1.2. 
The  discharge  curves  at  different  cycles  are  plotted  in  Fig.  13.  A  loss 
in  capacity  is  observed  immediately  after  the  first  cycle  indicating 
the  formation  of  SEI  film.  The  film  grows  with  time  increasing  the 
resistance  which  manifests  as  voltage  drop  with  increasing  cycle 
number.  The  increase  in  resistance  however  lowers  the  over¬ 
potential  for  side  reaction  which  slows  the  growth  of  the  film.  This 
leads  to  a  continuous  albeit  slower  loss  in  capacity  with  increasing 
cycle  numbers. 

The  SEI  film  thickness  as  a  function  of  time  is  presented 
in  Fig.  14.  The  film  growth  is  found  to  be  linearly  increasing  with  the 
cycling  time  and  hence,  in  extension,  it  is  linear  with  the  number  of 
cycles.  Brousselyet  al.  [41]  study  the  aging  of  lithium  ion  cells  under 
various  operation  regimes.  In  the  regime  where  the  conductivity  of 
the  SEI  film  becomes  rate  controlling,  they  propose  a  linear  scaling 
between  the  SEI  film  growth  and  the  number  of  cycles.  In  the 
present  work,  the  SEI  film  influences  the  cell  potential  through  an 
Ohmic  drop  across  it.  That  is,  the  overpotential  driving  the  elec¬ 
trochemical  reaction  is  directly  affected  by  the  film  resistance. 
Therefore,  the  scaling  law  provided  by  Broussely  et  al.  [41]  cor¬ 
roborates  the  linear  film  growth  observed  in  the  present  work.  This 
provides  a  consistency  check  for  the  description  of  degradation 
used  in  the  model. 

3.4.1.  Rate  dependence 

The  results  from  the  ROM  degradation  model  are  used  to  study 
the  rate  dependence  of  the  SEI  film  growth.  For  a  fixed  number  of 


Fig.  14.  The  growth  of  the  SEI  film  as  a  function  of  the  cycling  time  in  hours. 
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Fig.  15.  The  growth  of  the  SEI  film  varying  as  a  function  of  the  current  rate  for  a  fixed 
number  of  cycles. 

cycles,  the  increase  in  the  SEI  film  thickness  is  plotted  as  a  function 
of  the  current  rate  in  Fig.  15.  As  the  SEI  formation  occurs  during  the 
charge  cycle,  the  charge  rate  is  varied  while  keeping  the  discharge 
rate  constant  at  1C.  From  Fig.  15,  it  can  be  observed  that  for  a  fixed 
number  of  cycles  the  SEI  growth  decreases  with  increase  in  the 
charge  rate.  At  5C,  the  SEI  growth  at  the  end  of  5  cycles  is  almost 
equal  to  the  growth  observed  in  a  single  charge  cycle  at  0.1  C. 

This  rate  dependence  of  the  SEI  film  thickness  originates  from 
an  interplay  of  the  side  reaction  rate  and  the  time  available  for  its 
growth.  At  higher  current  rates,  a  constant  current  charge  cycle 
occurs  in  1  /5  th  of  the  time  it  would  take  at  1 C.  This  gives  lesser  time 
for  the  side  reaction  to  occur  and  hence  leads  to  lower  SEI  film 
growth. 

For  high  power  applications  low  cell  resistance  is  preferred.  This 
translates  to  minimal  SEI  growth  during  the  operation.  While  for 
high  energy  applications,  to  avoid  continuous  loss  of  lithium  to  the 
side  reactions,  formation  of  a  stable  SEI  in  the  first  cycle  followed  by 
a  quick  saturation  in  its  growth  is  preferable.  The  insights  from  the 
present  model  results  can  be  used  to  design  suitable  charging 
profiles.  For  example,  incorporating  a  higher  rate  constant  current 
charging  segment  can  lead  to  a  lower  SEI  growth  for  a  high  power 
cell.  While  a  lower  rate  constant  current  charging  can  lead  to  a 
rapidly  saturating  first  cycle  SEI  growth  for  a  high  energy  cell. 

3.5.  Operating  guidelines 

In  the  previous  sections,  the  implication  of  the  phase  change 
phenomenon  on  the  fresh  and  aged  cell  response  is  discussed.  The 
features  observed  as  a  result  of  the  two-phase  coexistence  are  used 
to  suggest  the  following  operating  guidelines  for  maximizing  uti¬ 
lization  of  cells  such  as  LFP/C 

•  For  a  fresh  cell,  utilization  can  be  maximized  if  similar  stages 
(charge  or  discharge)  succeed  one  another. 

•  To  obtain  a  lower  SEI  growth  (suitable  for  power  cells)  a  higher 
rate  constant  current  charging  segment  can  be  incorporated. 

•  For  a  saturated  SEI  film  (suitable  for  energy  cell)  a  lower  rate 
constant  current  charging  is  recommended  which  can  lead  to  a 
rapidly  saturating  first  cycle  SEI  growth. 

4.  Conclusion 

A  physics  based  reduced  order  model  (ROM)  is  proposed  for 
lithium  ion  cells  with  a  phase  change  electrode.  A  simple  and  novel 
approach  is  adopted  to  model  the  multi-phase  coexistence  within 
the  active  material  using  a  concentration  dependent  diffusivity.  The 


proposed  ROM  is  validated  against  experimental  data  from  litera¬ 
ture.  It  is  demonstrated  that  the  model  can  exhibit  the  character¬ 
istic  experimental  behavior  of  the  phase  change  electrodes.  It  is  also 
used  to  gain  physical  insight  into  various  features  of  a  phase  change 
system  such  as  charge  discharge  asymmetry  and  path  dependence. 
As  the  model  innately  takes  into  account  the  formation  of  multiple 
phases  it  enables  studying  the  cycle  life  of  cells  with  phase  change 
electrodes.  To  this  end,  the  effect  of  SEI  formation  due  to  side  re¬ 
actions  is  implemented  in  this  ROM  framework  and  validated  with 
experimental  data.  It  is  demonstrated  that  higher  rates  of  constant 
current  charging  lead  to  lower  film  formation.  Optimum  operating 
protocols  which  can  reduce  either  the  capacity  or  the  voltage  losses 
can  be  obtained  using  this  reduced  order  model.  To  conclude,  the 
present  ROM  based  on  physical  principles  can  be  used  for  BMS 
applications  due  to  its  algorithmic  structure.  It  can  also  be  used  to 
obtain  physical  insights  into  the  electrochemical  processes  leading 
to  optimum  utilization  of  the  cell.  Additional  physics  such  as 
thermal  effects  can  easily  be  adopted  into  this  framework  which 
forms  the  scope  of  future  work. 
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List  of  symbols 

a  specific  surface  area  of  the  porous  active  material,  m  1 

c  concentration  of  lithium,  mol  itT3 

cs  surface  concentration  of  lithium,  mol  nrT3 

D  diffusion  coefficient  of  lithium,  m2  s_1 

F  Faraday's  constant,  C  mol  1 

i  current  density,  A  irT2 

I  total  current  density,  A  m~2 

L  thickness  of  the  electrode/separator/cell,  m 

k  reaction  rate,  m  s_1 

j  reaction  rate  at  the  pore  wall  interface,  mol  m~2  s_1 

M  molecular  weight  of  the  SEI  film,  kg  mol-1 

q  mass  flux  at  the  electrode— electrolyte  interface, 

mol  m~2  s_1 

r  radial  coordinate,  m 

R  radius  of  the  active  particle  within  the  electrode,  m 

Rg  Ideal  gas  constant,  J  (mol  K)_1 

Rf  SEI  film  resistance,  Q  m2 

S  surface  area  of  the  electrode,  m2 

t+  transference  number  of  the  lithium  ions 

t  time,  s 

T  temperature,  I< 

U  open  circuit  voltage,  V 

V  cell  voltage,  V 

x  Cartesian  coordinate 

ac,  aa  charge  transfer  coefficients 
<5  SEI  film  thickness,  m 

e  porosity 

k  electrical  conductivity  of  the  electrolyte  phase,  S  m_1 

p  density  of  SEI  film,  kg  m~3 

a  electrical  conductivity  of  the  solid  phase,  S  nrT1 

4>  potential,  V 

Subscript 

f  film  variable 

in,  ip  interfacial  variables 

max  maximum  value  of  the  variable 

n  negative  electrode  properties 

ni,  pi  intercalation  reaction 
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p  positive  electrode  properties 

s  side  reaction  variable 

0  initial  value  of  the  parameter 

1  solid  phase  properties 

2  electrolyte  phase  properties 

a  lithium  deficient  phase 

(3  lithium  rich  phase 
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